function [fuzzy_center, fuzzy_member]=fuzzy_cmean(x,y0,c,epc,beta,N,isplot)
%introduction to parten recoginition: p/199
%it works fine for 2,4 dimension but nonsense for hingh dimension data??
%step 1
%initialization iteratioins
k=0;
%N=10; %the number of iterations
%beta=2; %the tuning parameter with controls the degree of fuzziness 
%epc=5*10^(-3); % the given tolerance
%c=2 % the number of c centers
%y0=[0 0; 2 2] % the initial c cluster centers is a nx1 vecator
%x =[ 0 0; 0 1; 1 1 ; 3 3 ;4 2]; % the given patterns 

y=y0;
[nrow ndim]=size(y0);
start=1;
if isplot==1
    %plot(x(:,1),x(:,2),'ro');
    clf;
    plot3(x(:,1),x(:,2),x(:,3),'bo');
    hold on;
    %plot(y(:,1),y(:,2),'bo');
    plot3(y(:,1),y(:,2),y(:,3),'ro');
    pause;
end
while start==1
    %step2 
    %calculate the distance
    d=sqrt(som_eucdist2(y,x));
    %step3
    %calcualte the membership
    sum_c_dist=sum(d.^(-2/(beta-1)),1);
    membership=[d.^(-2/(beta-1))./repmat(sum_c_dist,c,1)];
    %step4 updata centers
    %find di=0 and replace membership as 1 all other center's membership=0 
    [zeroD_i zeroD_j]=find(d==0);
    membership(:,zeroD_j)=0;
    lnZero=length(zeroD_i);
    for ki=1:lnZero
      membership(zeroD_i(ki),zeroD_j(ki))=1;
    end
    %updata senter have to be study more!
    sum_member=sum(membership.^beta,2);
    y2=membership.^beta*x./repmat(sum_member,1,ndim);
    %step5 check error
    fuzzy_error=0;
    for i=1:c
        fuzzy_error=fuzzy_error+som_eucdist2(y2(i,:),y(i,:));
    end
    fuzzy_error=sqrt(fuzzy_error);
    if fuzzy_error>epc 
        y=y2;
        it=k+1;
    else
        fuzzy_center=y;
        fuzzy_member=membership;
        break;
    end
    if k==N;
       disp('No convergence!  Need more N');
       fuzzy_error
       fuzzy_center=y;
       fuzzy_member=membership;
       break;
    else
        k=k+1;
        if isplot==1
            %plot(x(:,1),x(:,2),'ro');
            plot3(x(:,1),x(:,2),x(:,3),'bo');
            hold on;
            plot3(y(:,1),y(:,2),y(:,3),'ro');
            %plot(y(:,1),y(:,2),'bo');
        %y,d,membership
            pause;
        end
    end
end  %end all loop
